#include <arm_neon.h>
#include <stdlib.h>
#include <iostream>
#include <sys/time.h>
using namespace std;
//__attribute__((aligned(32))) float a[1001][1001];
void chuan(float** a, int n)
{
    int i, j, k;
    for (k = 0; k < n; k++)
    {
        for (j = k + 1; j < n; j++)
        {
            a[k][j] = a[k][j] / a[k][k];
        }
        a[k][k] = 1.0;
        for (i = k + 1; i < n; i++)
        {
            for (j = k + 1; j < n; j++)
            {
                a[i][j] = a[i][j] - a[i][k] * a[k][j];
            }
            a[i][k] = 0;
        }
    }
}
void bing_qi(float** a, int n)
{
    int i, j, k;
    for (k = 0; k < n; k++)
    {
        for (j = k + 1; j % 4 != 0 && j < n; j++)
        {
            a[k][j] = a[k][j] / a[k][k];
        }
        for (j; j <= n - 4; j += 4)
        {
            float32x4_t temp;
            float32x4_t t = vmovq_n_f32(a[k][k]);
            temp = vld1q_f32(a[k] + j);
            temp = vdivq_f32(temp, t);
            vst1q_f32(a[k] + j, temp);
        }
        for (j; j < n; j++)
        {
            a[k][j] = a[k][j] / a[k][k];
        }
        a[k][k] = 1.0;
        for (i = k + 1; i < n; i++)
        {
            for (j = k + 1; j % 4 != 0 && j < n; j++)
            {
                a[i][j] = a[i][j] - a[i][k] * a[k][j];
            }
            for (j; j <= n - 4; j += 4)
            {
                float32x4_t t1;
                float32x4_t t2;
                t1 = vld1q_f32(a[i] + j);
                t2 = vld1q_f32(a[k] + j);
                float32x4_t t3 = vmovq_n_f32(a[i][k]);
                t3 = vmulq_f32(t3, t2);
                t1 = vsubq_f32(t1, t3);
                vst1q_f32(a[i] + j, t1);
            }
            for (j; j < n; j++)
            {
                a[i][j] = a[i][j] - a[i][k] * a[k][j];
            }
            a[i][k] = 0;
        }
    }
}
void bing_bu(float** a, int n)
{
    int i, j, k;
    for (k = 0; k < n; k++)
    {
        for (j = k + 1; j < n; j += 4)
        {
            float32x4_t temp;
            float32x4_t t = vmovq_n_f32(a[k][k]);
            temp = vld1q_f32(a[k] + j);
            temp = vdivq_f32(temp, t);
            vst1q_f32(a[k] + j, temp);
        }
        a[k][k] = 1.0;
        for (i = k + 1; i < n; i++)
        {
            for (j = k + 1; j < n; j += 4)
            {
                float32x4_t t1;
                float32x4_t t2;
                t1 = vld1q_f32(a[i] + j);
                t2 = vld1q_f32(a[k] + j);
                float32x4_t t3 = vmovq_n_f32(a[i][k]);
                t3 = vmulq_f32(t3, t2);
                t1 = vsubq_f32(t1, t3);
                vst1q_f32(a[i] + j, t1);
            }
            a[i][k] = 0;
        }
    }
}
void m_rest(float** a, int n)
{
    for (int i = 0; i < n; i++)
    {
        for (int j = 0; j < i; j++)
            a[i][j] = 0;
        a[i][i] = 1.0;
        for (int j = i + 1; j < n; j++)
        {
            float x = rand();
            while (x == 0) x = rand();
            a[i][j] = x;
        }
    }
    for (int k = 0; k < n; k++)
        for (int i = k + 1; i < n; i++)
            for (int j = 0; j < i; j++)
                a[i][j] += a[k][j];
}
int main()
{
    srand((unsigned)time(NULL));
    int n;
    n = 10;

    while (n < 550)
    {
        float** a = (float**)aligned_alloc(32, n * sizeof(float**));
        for (int i = 0; i < n; i++)
            a[i] = (float*)aligned_alloc(32, n * sizeof(float*));
        m_rest(a, n);
        struct timeval tstart, tend;
        double timeUsed1 = 0.0, timeUsed2 = 0.0, timeUsed3 = 0.0;

        float** t1 = (float**)aligned_alloc(32, n * sizeof(float**));
        for (int i = 0; i < n; i++)
            t1[i] = (float*)aligned_alloc(32, n * sizeof(float*));

        gettimeofday(&tstart, NULL);
        chuan(t1, n);
        gettimeofday(&tend, NULL);
        timeUsed1 += 1000000 * (tend.tv_sec - tstart.tv_sec) + tend.tv_usec - tstart.tv_usec;


        float** t2 = (float**)aligned_alloc(32, n * sizeof(float**));
        for (int i = 0; i < n; i++)
            t2[i] = (float*)aligned_alloc(32, n * sizeof(float*));
        
        gettimeofday(&tstart, NULL);
        bing_qi(t2, n);
        gettimeofday(&tend, NULL);
        timeUsed2 += 1000000 * (tend.tv_sec - tstart.tv_sec) + tend.tv_usec - tstart.tv_usec;
        

        float** t3 = (float**)aligned_alloc(32, n * sizeof(float**));
        for (int i = 0; i < n; i++)
            t3[i] = (float*)aligned_alloc(32, n * sizeof(float*));

        gettimeofday(&tstart, NULL);
        bing_bu(t3, n);
        gettimeofday(&tend, NULL);
        timeUsed3 += 1000000 * (tend.tv_sec - tstart.tv_sec) + tend.tv_usec - tstart.tv_usec;

        cout << n << " " << (timeUsed1 / timeUsed2) * 100 << "%" << " " << (timeUsed1 / timeUsed3) * 100 << "%" << endl;
        if (n >= 200) n += 50;
        if (n < 200) n += 10;
    }
    return 0;
}